Anisotropy and universality: the Oslo model, the rice pile experiment and the 

quenched Edwards- Wilkinson equation. 
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We show that any amount of anisotropy moves the Oslo model to another known universahty 
class, the exponents of which can be derived exactly. This amounts to an exact solution of the 
quenched Edwards- Wilkinson equation with a drift term. We argue that anisotropy is likely to be 
experimentally relevant and may explain why consistent exponents have not been extracted in the 
rice pile experiments. 
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The suggestion that power-law like distributions fre- 
quently reported in experiments may be the effect of 
threshold dynamics and metastabihty was made by Bak, 
Tang and Wiesenfeld in 1987 Q and was called Self- 
Organized Criticality (SOC). Avalanche dynamics in 
granular piles was from the onset used as a metaphor 
and laboratory for SOC behavior and has inspired many 
models and experiments |2| . One of the most celebrated 
of these efforts is the experimental study of avalanches 
in one dimensional ricepiles by Frette et a/.0] and the 
theoretical Oslo model |j| inspired by the rice pile exper- 
iment. The general interest and relevance of such studies 
rely on the assumption, guided by equilibrium critical 
phenomena, that the critical behavior of scale invariant 
systems falls into universality classes determined solely 
by a few general characteristics of the system, such as 
symmetry and dimension. So-called "relevant" param- 
eters can decide which of the symmetries the system is 
asymptotically dominated by. 

The role of anisotropy in SOC has been highlighted 
very early by Hwa and Kardar |^ and Grinstein et al. , 
who used anisotropic Langevin equations to describe 
sandpiles. On the cellular automata level, Kadanoff 
et al. 7] have conjectured that the net flux of particles 
is a relevant parameter. In this paper we confirm this 
conjecture for the Oslo model, which shows a clear cut 
and consistent relevant dependence on anisotropy. This is 
of great importance for the interpretation of experimen- 
tal results Q and more generally for the much studied 
quenched Edwards- Wilkinson equation |^ 0] . 

Moreover, contrary to suggestions in former stud- 
ies, the switch between different universality classes 
(crossover) is not t rigg ered by the introduction of 
stochasticit y Eol llll ll2L ll^ nor by multiple topplings 
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Similar to |l7l |. the system size at which the crossover 
occurs, Lx, depends on the strength of the anisotropy 
V. We exemplify two possible mechanisms causing 
anisotropy in experiments, one of which vanishes with 
system size L fast enough to keep vL constant. This rep- 
resents a marginal case and consequently makes a unique 
identification of the critical exponents impossible. 



The original Oslo model (oOM) consists, in one di- 
mension, of a lattice of sites i = 1, . . . , L. Two coupled 
dynamical variables are associated with each lattice site: 
the primary variable Zi G {0, 1,2,.. .} and the threshold 
variable zf £ {1,2}. The initial configuration consists of 
Zi = \fi and a random configuration of the zf. The 
system is driven by increasing zi by one (a "grain" ) fol- 
lowed by a relaxation of all sites 1 < i < L for which 
Zi > ("over-critical" sites). In case a site i is over- 
critical, the following updates are performed ("toppling" 
or "relaxation"): Zi —>■ Zi — 2 and Zi±i — > Zi±i + 1 and, 
importantly, the existing value of the threshold zf is af- 
terwards replaced by 1 with probability p and by 2 with 
probability 1 — p. The boundaries are updated the same 
way except that for i = 1 (i = L) addition on site 
(i -I- 1) is omitted. 

We now introduce a tunable degree of anisotropy into 
the dynamics. An over-critical site i is relaxed in the 
following way. Only left movement: with probability 
pi{l ~ Pr) perform the updates: Zi ^ Zi — 1 and Zi-i 
Zi-i + 1. Only right movement: with probability Pi) 
perform the updates Zi Zi — 1 and z^+i —>■ z^+i -I- 1. 
Both left and right movement: with probability pip^ per- 
form the updates Zi±i Zi±i + 1 and Zi ^ Zi — l. A new 
zf is chosen, at random, after every successful update, 
i.e. when at least one grain has been redistributed. We 
call this version of the model the anisotropic Oslo model 
(aOM). The strength of the anisotropy is described by the 
drift velocity v = (pr —pi)/ {pr + Pi ) which is the net flux 
of grains through the system. Clearly it is only sensible to 
study Pr + pi > 0. The case pr = pi = 1 corresponds ex- 
actly to the oOM, while p,. — pi ^ 1 represents a stochas- 
tic variant of the oOM. The avalanche exponents for the 
extreme, totally asymmetric case pi = and p,. = 1 can 
be obtained exactly 0| and describe, as we shall see be- 
low, the scaling behavior for all v > 0. We are interested 
in the statistics of the sizes, s, of the avalanches of re- 
laxation induced by the driving zi zi + 1. The size 
of an avalanche is in both versions of the model deflned 
as the number of times the relaxation rule was success- 
fully applied after the drive zi — > zi + 1 in order to make 
Zi < zf Vi yet again. Thus s > 0. From the definition 
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above it is clear that the model is equally well defined in 
terms of reduced probabilities, omitting the case where 
no grain is redistributed. 

In the totally anisotropic or asymmetric limit |0| the 
model resembles some features of other exactly solved, 
directed models 0, 0, 0, 0, 0| . In two dimensions 
a very similar model has been studied numerically [l4l| . 
However, we stress that contrary to some other "exact" 
solutions, the model is solvable directly on the lattice 
and without assuming any scaling behavior 'is] . Also the 
amplitudes of the moments can be calculated exactly. 

It is a tedious, but straightforward task to show that 
the aOM is "Abelian", i.e. the order of updates is ir- 
relevant for its statistical properties. Since the micro- 
dynamics which prescribes the order of updates is irrel- 
evant, there is no unique way to define a microscopic 
timescale. Presumably universal exponents of the du- 
ration of avalanches are therefore mainly a property of 
the arbitrary choice of the microdynamics. According 
to Hughes and Paczuski a non- Abelian variant of 
an Abelian model may or may not remain in the same 
universality class. 

Moreover, one notes that the oOM as well as the aOM 
contains multiple topplings, i.e. a single site can relax 
several times during a single avalanche. 

We now describe the avalanche statistics of the aOM. 
In Fig. n we demonstrate that the avalanche size distri- 
bution V{s] L,pr,pi) follows simple (finite size) scaling 



P{s;L,pr,pi) = aipr,pi)s '^Q 



b{Pr,Pi)LD 



for s > si 
(1) 

where Q is the universal scaling function, a{pr,pi) and 
b{Pr,Pi) are two anisotropy and system dependent pa- 
rameters, and si is the lower cutoff independent of L. 
The values for the scaling exponents are Ta = 4/3 and 
Da = 3/2, which can be derived exactly in the asymmet- 
ric limit t; = 1 and represent the known universality 
class of the stochastic, directed sandpile in two dimen- 
sions 0, 0, 0| , which is in turn closely related to the 
directed sandpile . Numerically, these exponents have 
been found for all f > studied at sufficiently large sys- 
tem sizes L ^ Lx- The crossover that occurs around 
Lx is discussed in detail below. The two exponents of 
the aOM are to be compared with the exponents for the 
oOM, of To = 1.556(4) and = 2.25(2). Since the av- 
erage avalanche size scales linearly with the system size, 
the exponents are related by D{2 — r) = 1. \4. l21| 

The easiest way to derive the exponents from numeri- 
cal data is by analysis of the moments [l^l, which scale 
according to 1^ for n > t — 1 in leading order like 

/•OO 

(s") = / dss^'Vis; L,pr,pi) = a(6L^)i+"-^g„ + . . . 
Jo 

(2) 

where (?„ is discussed below and . . . denotes subleading 
terms, especially Wegner's corrections to scaling 23]. In 




FIG. 1: Main panel: Data collapse {Vis; L,pr,pi)s'^ vs. 
s/L^) of the normalized and binned Vis; L,pr,pi) for the 
two extreme cases (p,- = 1 with pi — and pi = 1) for 
L = 1280,2560,5120. The PDF r{s; L,pr,pi) has been 
rescaled using r — 1.333. . . and D = 1.5 for the aOM and 
T = 1.555 . . . and D = 2.25 for the oOM. Inset: Rescaled 
PDFs for three different choices of pr > Pi, namely {pr,pi) = 
(1.0, 0.95), (1.0, 0.25), (0.75, 0.25), and L = 640,1280,2560, 
and two choices of pr ~ Pi, namely 0.75 and 0.25 and the 
same range of L. The rescaled PDFs of each tuple {pr,Pi) 
would form a single line, which fixes the two exponents r and 
D. By tuning the parameters a and b in Eq. Q, as done in 
the inset, the resulting single line collapses with one of the 
extreme cases plotted in the main panel and the inset. De- 
viations for intermediate values of s/L^ (i.e. away from the 
bump towards smaller values of s) are expected to vanish at 
sufficiently large L. 
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FIG. 2: Scaling of (s^)/L^'^ for different anisotropies. The 
filled circles are results for pr — pi as indicated and pr = 
1. Open circles show other parameters {pr,pi)- The dashed 
lines are the two extreme cases oOM {pr = pi = 1) and the 
solvable model {pr = l,pi = 0). The open arrows mark the 
approximate crossover. 
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FIG. 3: Scaling of (see Eq. (|3J) for different anisotropies. 
The filled circles are results for pr — pi in the same vertical 
order as shown in Fig.|5|and pr = 1. Open circles show other 
parameters {pr,Pi) as indicated. The dashed lines are the 
two extreme cases oOM {pr — pi = 1) and the solvable model 
(pr = ~ 0). The open arrows mark the crossover points 
in Fig. 11 



the following, the crossover is studied by means of the 
rescaled second moment, (s^)/L^/^, which is shown in 
Fig. 121 For non- vanishing anisotropy, w > 0, it ap- 
proaches a constant as L —^ oo. For very small but 
finite values of v and L the rescaled moment increases 
with L like L" ''^, corresponding to the oOM behavior, 
but at L « Lxiv) it crosses over and eventually con- 
verges to a finite constant. Below, we shall relate the 
behavior of Lx{v) to the effective anisotropy relevant to 
an experiment of a given size. Here we emphasize that 
Fig. 13 clearly demonstrates that the universality class of 
the extremely anisotropic case, Pr = 1, pi = 0, contains 
all systems with non-vanishing anisotropy v > 0. That 
renders the oOM with v = a special case; however, it 
is remarkable that even for pr = pi ^ 1 the model still 
shows oOM behavior. Thus, it is not the stochasticity 
itself [lol im which induces the change in critical 
behavior. 

Eq. 1^ allows the definition of universal amplitude ra- 
tios 



is') 



(n-l) 



(3) 



which can easily be proven to be asymptotically inde- 
pendent of a, b and L. The two constraints on Q, which 
fix a and 6 in can be chosen such that g„ are the 
moments of x~'^Gix) as used in namely by imposing 
that J^x'^-'^Gix) = J^x^-^g{x) = 1. The universal 
amplitude ratio gs as shown in Fig. |3 indicates not only 
the same crossover behavior as observed in Fig. |21 but 
also the universality of Q. 



The importance of the above result is highlighted when 
we recall that the oOM in the continuum limit is de- 
scribed by the quenched Edwards- Wilkinson (EW) equa- 
tion 0. A similar derivation shows that the aOM is a 
quenched EW equation with an additional drift term 

dth{x, t) = Ddlh[x, t) - vd^h{x, t) + r]{x, h{x, t)), (4) 

where D is the surface tension and v the anisotropy or 
drift velocity as defined above. In the aOM, the quenched 
noise r]{x, h{x, t)) represents the randomly chosen zf and 
h{x, t) is the number of charges received by site x at time 
t. The quenched nature of the noise makes it difficult to 
solve (0J) directly. Together with the boundary conditions 
0, it prevents the drift term vd^h from being absorbed 
by a Galilean transformation. However, the above results 
determines the roughness exponent via D = 1 + x to 
be X = 1/2 for u > 0, which has already been suggested 
in another case of anisotropic depinning j3| ■ For v = 
numerics for the oOM suggest correspondingly that x = 
1.25(2). 

The crossover for small v > with increasing L can 
also be illuminated by a study of the individual grains 
in the system, which behave like biased random walk- 
ers. This leads again to a diffusion equation with drift 
term and two absorbing boundaries ,25j. The average 
avalanche size is the average time the particles spend in 
the system divided by the average number of grains re- 
distributed per toppling. The crossover is expected as 
soon as the ballistic motion dominates over the diffusion, 
L'^/D > L/v, thus Lx{v) = D/v. This has also been 
tested numerically, based on heuristic estimation of Lx, 
as shown by the marks in Fig. |21 
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FIG. 4; A stylized toppling of a single grain: If an elongated 
grain of width b and height a topples from site i to site i + 
it reduces the height at i by fe and increases the height at i -f 1 
by a, thereby increasing the slope at site i — 1 by 6, decreasing 
the slope at i by a -I- 6 and increasing the slope at i -I- 1 by a. 
The net flux of slope is a — b to the right. 

We now discuss the relevance of anisotropy to real 
experimental granular systems. We stress that the 
anisotropy is in the amount of slope transported between 
sites involved in a relaxation event. Any net flux of the 
slope is eventually compensated by the toppling of the 
last site and the slope is therefore asymptotically sta- 
tionary. 
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One process leading to an anisotropic redistribution 
of slope arises when the toppling grain is elongated, see 
Fig. ^ Most remarkably, in the original experiment Q 
it was noted that only the elongated rice samples showed 
scale invariant behavior. A reorientation of a single grain, 
as shown in Fig.^ leads to a net flux of slope to the right. 
It can happen only once, and in fact it depends on how 
and whether the rice enters and/or leaves the system with 
a typical orientation. If that is the case, then an average 
reorientation is distributed among all topplings on its way 
through the system, i.e. v oc 1/L. Because Lx oc 1/v 
and V oc 1/L, this represents a marginal case, because it 
is impossible to decide whether L Lx{v) or not. The 
dependence of, say, the ratio (s^)/L^-^ on L would in 
this case be given by an (unknown) trajectory through 
the diagram in Fig. |21 since a change in L also leads to 
a change in v, which allows non- universal quantities to 
enter. 

We do not know if the experiment by Frette et aL0 
involves this complication, but the exponents ex- 
tracted from the experiment by the authors are not 
consistent . 

Local rearrangements, such as expansions (on the site 
that looses a grain) and compressions (on the site that 
receives the grain) lead to an anisotropy which does 
not vanish with L: Say column i — 1 from which the 
grain leaves expands by e{hi) and column i + 1 which 
receives the grain is compressed by an amount e(ft.i+i). 
Then the changes of the slopes during toppling are 
Azi_i 1 - e{hi), Azj = -2 + e{hi) + e{h,+i) and 
Azi+i = 1 — e(/ii+i). Assuming that the columns behave 
elastically, e would be an increasing function of h which 
results in a net flow to the right, i.e. v > 0. However, it 
remains unclear whether any of these effects can be seen 
in experimental systems. Moreover, having shown that 
anisotropy is a relevant field, it would not be surprising 
to find other relevant fields which lead to yet another 
universality class. 

We have demonstrated that for any amount of 
anisotropy the exponents of the original Oslo model 
change and are given by simple rational numbers which 
can all be obtained exactly [ig- The crossover has been 
studied numerically using a moment analysis, Eq. 
and universal amplitude ratios, Eq. and the crossover 
length has been determined. The generalized model de- 
scribed above continuously connects the established orig- 
inal Oslo model and an exactly solvable, directed variant. 
This variant has, compared to the original Oslo model, 
an enormous basin of attraction, so that the latter may 
be regarded a special case of the former. 

Moreover, we find a change in critical behavior of an 
SOC model, genuinely due to anisotropy, rather than 
stochasticity or the presence of multiple topplings. 

The results are theoretically interesting especially be- 
cause of their relation to the Edwards- Wilkinson equa- 
tion, the roughness exponent of which has been obtained 



in case of the presence of a drift term in one dimension 
to be X = 1/2. Moreover, according to our study, exper- 
iments are seriously complicated due to a coupling be- 
tween system size and effective anisotropy. That might 
provide a clue as to the apparent difficulties to find the- 
oretically predicted exponents in the real world. 
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